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ON THE STABILITY OF NON FORCE-FREE MAGNETIC EQUILIBRIA IN STARS 
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ABSTRACT 

The existence of stable magnetic configurations in white dwarfs, neutron stars and various non- 
convective stellar regions is now well recognized. It has recently been shown numerically that vari- 
ous families of equilibria, including axisymmetric mixed poloidal-toroidal configurations, are stable. 
Here we test the stability of an analytically-derived non force-free magnetic equilibrium, using three- 
dimensional magnetohydrodynamic simulations: the mixed configuration is compared with the dy- 
namical evolution of its purely poloidal and purely toroidal components, both known to be unstable. 
The mixed equilibrium shows no sign of instability under white noise perturbations. This configu- 
ration therefore provides a good description of magnetic equilibrium topology inside non-convective 
stellar objects and will be useful to initialize magneto-rotational transport in stellar evolution codes. 
Subject headings: magnetohydrodynamics (MHD) — stars: interiors — stars: magnetic field — stars: 
neutron — Sun: magnetic topology — white dwarfs 



1. INTRODUCTION 

Magnetic fields are being detected more and more rou- 
tinely at the surface of many stars, and are responsible 
for various physical phenomena likely to deeply modify 
our traditional vision of stellar evolution, especially 
during their early and late stages. The presence of 
even a relatively weak magnetic field can have an 
important impact on the collapse and f ragmentation 
of prestellar cores (jCommercon et al1l2010f ). as well as 
influenci ng the rotation rate of pre-main-sequence stars 
(see e.g. lAlecian et al.l l2008h . On the other side of the 
Hertzsprung-Russell diagram its feedback effects may 
also play a key role in supernovae and mechanical energy 
deposition in the interstellar medium, for instance. 
Magnetic fields are also an important actor in main- 
sequence stars. First, we cannot dismiss the possibility 
of a large-scale magnetic field being responsible for the 
quasi-uniform rotation behaviour in the bulk of the solar 
radiation zone, as reveal ed by p- modes helioseismology 
(|Eff-Darwich et al.l l2Q08h . Second, strong fields (300 
G to 30 kG) are observed via the Zeeman effect in 
som e fraction of main-s equence A stars (the Ap stars, 
see lAuriere et al.l l2007[ h as well as in B sta rs and 
in a handful of O stars (jGrunhut et al.l l2009f ). The 
bimodality of rotational velocities observed among Ap 
vs. normal A stars shows the critical effect of magnetic 
fields on rotation and therefore also on meridional 
circu lation and chemical transport (see iMathis fc Zahnl 
2005). Finally, magnetic white dwarfs display fields 
strength of 10 4 — 10 9 G, and neutron stars host fields 
in the range 10 8 — 10 15 G, in both cases detected using 
several distinct methods. 

The large-scale, ordered nature of these fields (often 
approximately dipolar) and the scaling of their strengths 
as a function of their host properties (according to the 
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flux conservation scenario) favour a fossil hypothesis, 
whose origin is not yet elucidated. 

Another fundamental question is the topology of 
these large-scale magnetic fields. To have survived 
since the star's formation, a field must be stable on 
a dynamic ( Alfve n) timescale. It was suggested by 
iPrendergastl (|1956f ) that a stellar magnetic field in sta- 
ble axisymmetric equilibrium must contain both poloidal 
(meridional) and toroidal (azimuth al) components, since 
both are unstable on their own (lTavlerl)l975 [Wrightl 
119731 : iBraithwaitd 120061 : iBonanno fc Urpinll2008ri ). This 
was confirmed recen tl y by n umerical simulations by 
Brait hwaite fc Spruitl ([20041 ): iBraithwaite fc Nordlundl 
(2006|) who showed that an arbitrary initial field evolves 
on an Alfven timescale into a stable configuration; ax- 
isymmetric mixed poloidal-toroidal fields were found. 
Once formed, it continues to evolve on longer timescales 
through diffusive processes such as finite conductivity: 
the field then moves outwards, passing quasi- statically 
through a series of stable axisymmetric equilibria until 
it changes eventually to a non- axisymmetric equilibrium. 
These n on- axisymme t ric e quilibria are described in more 
detail in IBraithwaitd |2009). 

2. THE RELAXED NON FORCE-FREE CONFIGURATION 

Here, we deal with axisymmetric, non force- free mag- 
netic configurations (i.e. with a non-zero Lorentz force) 
in equilibrium inside a conductive fluid. We first restrict 
ourselves to the non-rotating case, but results also ap- 
ply t o rotating stars where rotation is uniform ([Woltjerl 
1959), which could be the case if magnetic field is strong 
enough, and where meridional circulation can be ne- 
glected (i.e. when the star does not loose angular mo- 
mentum and have a stationary st ructure: see[Busse 1981; 
Zahn 19921: iDecressin et al.ll2009l ). The more general case 
including (differential) rotation (and induced meridional 
circulation) will be studied in a forthcoming work. Sev- 
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eral reasons inclined us to focus on non force-free equi- 
libria instead of fo rce- free ones ; let u s briefly describe 
them here. First, iReisenegger (2009) reminds us that 
no configuration can be force- free everywhere. Although 
there do exist "force- free" configurations, they must be 
confined by some region or boundary layer with non-zero 
or singular Lorentz force. Discontinuities such as current 
sheets are unlikely to appear in nature except in a tran- 
sient manner. Second, non force-free equilibria have been 
identified in plasma physics as the result of relaxation 
(self-organization process invol ving magnetic reconnec- 
t ions, in resistive MHD ), e.g. by iMontgomerv fc Ph illips 
(11988ft: iShaikh et al.l ([2008ft . Third, as shown by 
iDuez fc Mathisl ((2010), this family of equilibria is a gen- 
eralization of Taylor states (force-free relaxed equilib- 
ria) in a stellar context, where the stratification of the 
medium plays a crucial role. 

2.1. The magnetic field in MHS equilibrium 

Let us briefly recall the assumptions made in building 
the semi-analytical mode l of magnetoh ydrostatic (MHS) 
equilibrium described bv lDuez fc Mathisl (|20 10ft . The ax- 
isymmetric magnetic field B(r, 9) is expressed as a func- 
tion of a poloidal flux \£(r, 0), a toroidal potential F(r, #), 
and the potential vector A (r, 0) so that it is divergence- 
free by construction: 



B = 



1 



r sin t 



(V* xe (p +Fe (p ) = V x A, 



(1) 



where in spherical coordinates the poloidal component 
is in the meridional plane (e r ,e<9) and the toroidal com- 
ponent is along the azimuthal direction (e^). The MHS 
equation expressing balance between the pressure gradi- 
ent force, gravity and the Lorentz force is 

= - V P - p VV + — (V x B) x B, (2) 

Po 

where V is the gravitational potential which satisfies the 
Poisson equation : V 2 V = 4irGp. 

2.2. The non force- free barotropic equilibrium family 
2.2.1. A variational approach 

Here, we focus on the minimum energy non force- free 
MHS equilibrium in stably stratified radiation zones. 
Given the field strengths in real stars, the ratio of the 
Lorentz force to gravity is very low: stellar interiors are 
in a regime where /3 = P/Pm^ » 1, ^Mag = B 2 / (2/io) 
being the magnetic pressure. We then identify the two 
MHD invariants governing the evolution of the recon- 
nection phase, that leads to relaxed states in the non 
force- free context: the magnetic helicity H = j v A-BdV 
and the mass encompassed in poloidal magnetic surfaces 
= J v ^f ~p dV, conserved because of the stable strat- 
ification. Assuming a selective decay during relaxation 
(the magnetic energy decays much faster than % and My 
so that they can be considered constant on an energetic 
decay e- folding time), a variational method allows us to 
derive the elliptic linear partial differential equ ation gov- 
erning (jWoltierlll959t buez fc Matriis1[20Toft : 
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Here, p is the density in the non-magnetic case, A*\I> = 
d rr ^ + sin# do {8q^ / s'mO) /r 2 the Grad-Shafranov op- 
erator in spherical coordinates, Ai a coefficient to be 
determined, R a characteristic radius, and /3q is con- 
strained by the field's intensity. This equation is similar 
to the Grad-Shafranov equation u sed to find equilibria 
in magnetically confined plasmas (jGrad fc Rubinl ll958: 
Shafranov 1966), the source term being here related to 
the s tellar structure through ~p (see iHeinemann fc Olbertl 
I1978L for a discussion of the general form of this equa- 
tion in astrophysics). Furthermore, this equilibrium is 
in a barotropic state (in the hydrodynamic meaning of 
the term, i.e. isobar and iso-density surfaces coincide) 
where the field is explicitly coupled with stellar struc- 
ture through: V x (Fc/p) = 0, where Fc is the Lorentz 
force. This is a generalization of Prendergast's equilib- 
rium taking into ac count compres sibility, first studied in 
polytropic cases bv lWoltierl (fl960ft . 

2.2.2. Solution 

The boundary conditions have now to b e discussed. In 
IDuez et al.l (|2010ft : lDuez fc Mathisl (|2010ft . we considered 
the general case of a field confined between two radii, 
owing to the possible presence of both a convective core 
and a convective envelope and to ensure the conserva- 
tion of magnetic helicity. We here choose to cancel both 
radial and latitudinal fields at the surface, to avoid any 
current sheets, conserving once again magnetic helicity. 
Owing to its small extension, the possible effects of the 
convective core on the large-scale surrounding field are 
neglected. Using Green's function method we finally ob- 
tain the purely dipolar, general solutions indexed by i: 
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R being the upper boundary confining the magnetic field; 
A^ are the set of eigenvalues indexed by i allowing to ver- 
ify the boundary conditions. The functions ji and y\ are 
respectively the spherical Bess el functions of the firs t and 
the second kind. As shown in lDuez fc M athis (2010]), the 
first radial mode is the lowest energy state. We thus focus 
here only on this mode i = 1. The toroidal magnetic field 
is then given using F($?) = Xi^f/R; in the case of a stably 
stratified n = 3 polytrope and for the simulation purposes 
where we set R = 0.85 i?*, we have Ai ~ 32.95, while for 
a constant density profile (in a zero gravity medium), we 
have Ai ~ 5.76. The solution for the n=3 polytrope is 
represented in Fig. [H The ratio of the poloidal to 
total magnetic energy density B?/B 2 (r, 6) is plotted as 
a function of the radius in Fig. for various latitudes. 
Notice that its integrated value E p /E ~ 5.23 x 10 -2 is in 
the range of ratios found in stable axisymmetric equilib- 
ria forming in simulations from random sm all-scale initial 
conditions (jBraithwaite fc Nordlu nd 2006). 

3. STABILITY: NUMERICAL METHOD 

3.1. The numerical model 

The setup of the numerical model is similar to that 
in Braithwaite & Nordlund 2006, where a fuller account 
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Fig. 1. — Toroidal magnetic field strength in colorscale (arbitrary 
field's strength) and normalized isocontours of the poloidal flux 
function ^ in meridional cut for the first equilibrium configuration 
(Aj; ~ 33). The neutral line is located at r c=l 0.23 R*. 
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Fig. 2. — Poloidal to total magnetic energy density B 2 /B 2 as a 
function of the radius for different colatitudes, and radial profile of 
the poloidal flux function ^ at the equator. 

can be found; a brief outline is given here. W e use the 
Stagger code ([Nordlund fc Galsgaardl IT995). a high- 
order finite-difference Cartesian MHD code containing 
a "hyper-diffusion" scheme. We use a resolution of 192 3 . 

We model the star as a self-gravitating ball of ideal 
gas (7 = 5/3) with radial density and pressure profiles 
initially obeying the polytropic (thus barotropic) relation 
P oc p 1+ ( 1//n \ with index n = 3 - a good approximation 
to an upper-main-sequence star. It seems unlikely that 
a different EOS, for instance that of a neutron star, will 
make even much quantitative difference to the results. 
The important point is the stable stratification. 

We use this model to compare the dynamical evolution 
of the mixed (poloidal-toroidal) configuration to that of 
its purely poloidal and toroidal components on their own, 
both of which are unstable as mentioned above. We 
should therefore see these instabilities, growing on an 
Alfven timescale. To test the stability of the configura- 
tions, we add a random "white noise" perturbation to the 



density field. The perturbation in density (1% in ampli- 
tude) contains length scales ranging from R* to 0.08i?*, 
the latter being double the Nyquist wavelength. This 
is roughly equivalent to azimuthal wavenumbers up to 
m = 38 at a radius of R*/2. 

3.2. Results 

Purely poloidal component — The simulation is run for 
around ten Alfven crossing times ta, over which time 
the instability grows, becomes nonlinear and results in 
the destruction of most of the original magnetic energy. 
The magnetic field amplitude is plotted at the top-left 
of Fig. [3l split into components according to azimuthal 
wavenumber m; obviously at t = all the energy is in the 
axisymmetric m = part. The bottom- left plate of the 
figure shows the mean velocity in each azimuthal mode. 
Note the clear transition at t « 2ta from the linear phase 
to the nonlinear, reconnective phase. 

Purely toroidal component — The middle plates of FlG. [3] 
show the evolution of the toroidal field - clearly, the 
m = 1 mode is dominant. In the lower panel of Fig. [4] 
are drawn the magnetic field lines for the evolved con- 
figuration, i.e. at t « 10 ta- After the linear growth of 
the m = 1 mode, the Tayl er instability manifests itself in 
the nonlinear regime (cf. [Bnml [20071 : lElstner et afl feOOS) 
mainly in the movement of spherical shells relative to 
one another - which is simulated here for the first time. 
We expect eventual complete destruction of the field on 
a longer timescale; however a more detailed investigation 
is left to a forthcoming paper. 

Mixed configuration — The mixed poloidal-toroidal 
configuration exhibits completely different behaviour. 
The magnetic and velocity amplitudes are plotted on 
the right of Fig. O where we see an absence of growing 
modes. The kinetic energy present results simply from 
the initial perturbation and the oscillations and waves it 
sets up. In Fig. [4] are drawn the magnetic field lines at 
t = and 10 ta; no significant change is seen. 

To better examine the poten tially unstab le regions, 
we use Tayler's stability criteria ([Tavlerll 19731 ) for purely 
toroidal fields and estimate th e stabilisation from the 
poloidal component, following iBraithwaitd (|2009). In 
Fig. [5]we plot Tayler's criteria for modes m = and 
rn = 1 the rn = mode is unstable almost everywhere 
and the m = 1 mode is unstable in a large cone 
around the poles; however the poloidal field stabilises 
these modes in most of the meridional plane except 
near the equatorial plane where it merely stabilises 
all wavelengths small enough to fit into the available 
space. We can examine closely the behaviour of the field 
in the vicinity of the magnetic axis, where it can be 
approximated as the additio n of an axial and a t oroida l 
field (cylindrical geometry). Bo nanno fc Urpinl (|2008a) 
outlined that in this case magnetic configurations can be 
subject to non- axisymmetric resonant instability. They 
determined the dependency of the Tayler instability 
maximum growth rate as a function of the azimuthal 
wave-number m and of the ratio e of the axial field to 
the toroidal one. In our case, close to the center the 
flux function exhibits a behaviour in ^ ex r 2 , so the az- 
imuthal field is proportional to s = r sin 6 corresponding 




Fig. 3. — Time evolution of the (log) amplitudes in azimuthal modes m = to 4 averaged over the stellar volume of the magnetic field 
(top row) and the velocity field (bottom row) in the simulations with the purely poloidal field (left), purely toroidal field (middle) and the 
mixed field (right). Initially, all the magnetic energy is in the m = mode since the initial conditions are axisymmetric. 



Fig. 4. — Magnetic field lines 
representing the mixed field con- 
figuration (right), looking along 
the axis (top) and from the 
side (middle and bottom). The 
purely toroidal component of the 
field is represented on the left. 
The colorscale is a function of 
the density. Upper and middle 
panels : configurations at t = 0; 
Lower panels : configuration at 
t = 10r A . 





to the Bonanno et al.'s parameter a = 1. As underlined 
by the authors, in that case the maximum growth rate 
changes remarkably slowly with m for all modes with 
m ^ 2 and the instability is weakly non-anisotropic. If 
we take as a value for s\ the radius of the neutral line or 
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Fig. 5. — Half of the meridional plane, showing the regions stable 
against the m = and 1 Tayler modes in the absence of the poloidal 
component, and their stabilisation by the radial component B r . 

the one where the azimuthal field is strongest, we obtain 
respectivel y e = 0.64 or e = 0.79. According to their 
study (see iBonanno fc Urpinl l2008al Fig. 7), we fulfill 
the stability criterion for the modes m = 0, 1 and 2. 
Our results are therefore in agreement with their linear 
analysis. 

In the simulations we run, the mixed configuration 
has a poloidal energy fraction E p /E = 0.052. The 
magnetic-to-thermal energy ratio E/U ~ 1/400, which 
should mean that for stability we require E p /E > 0.04 
(|Braithwaitell2009l ). We see then that this value of E/U 
is near the upper limit for stability - in other words, 
we are near the boundary of validity of the weak-field 



5 



approximation used in Section [2j 

4. CONCLUSION 

Using semi-analytic methods we derived (with an ap- 
propriate choice of boundary conditions) then tested an 
axisymmetric non force-free magnetostatic equilibrium 
which could exist in any non-convective stellar region: 
the radiative core of solar- type stars, the external enve- 
lope of massive stars, and compact objects. Using numer- 
ical simulations, we demonstrate the ability of the set-up 
to recover well-known instabilities in purely poloidal and 
toroidal cases, then find stability of the mixed config- 
uration under all imaginable perturbations. We show 
the agreement of the result with linear analysis (lim- 
ited in perturbations), highlighting the stabilizing in- 
fluence of the poloidal field on the toroidal one, espe- 
cially in the region close to the symmetry axis where 
purely toroidal fields usually develop kink-type instabil- 
ities in priority. This is the first time the stability of 
an analytically-derived stellar magnetic equilibrium has 
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thermore it will help to appreciate the internal magnetic 
structure of neutron stars, and various astrophysical pro- 
cesses involving magnetars (intense activity in the X-ray 
and gamma-ray spectra, quasi-periodic oscillations and 
eventually gamma-ray bursts). 
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